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We consider the flavor evolution of a dense neutrino gas by taking into account both vacuum 
oscillations and self interactions of neutrinos. We examine the system from a many-body perspective 
as well as from the point of view of an effective one-body description formulated in terms of the 
neutrino polarization vectors. We show that, in the single angle approximation, both the many- 
body picture and the effective one-particle picture possess several constants of motion. We write 
down these constants of motion explicitly in terms of the neutrino isospin operators for the many- 
body case and in terms of the polarization vectors for the effective one-body case. The existence of 
these constants of motion is a direct consequence of the fact that the collective neutrino oscillation 
Hamiltonian belongs to the class of Gaudin Hamiltonians. This class of Hamiltonians also includes 
the (reduced) BCS pairing Hamiltonian describing superconductivity. We point out the similarity 
between the collective neutrino oscillation Hamiltonian and the BCS pairing Hamiltonian. The 
constants of motion manifest the exact solvability of the system. Borrowing the well established 
techniques of calculating the exact BCS spectrum, we present exact eigenstates and eigenvalues of 
both the many-body and the effective one-particle Hamiltonians describing the collective neutrino 
oscillations. For the effective one-body case, we show that spectral splits of neutrinos can be 
understood in terms of the adiabatic evolution of some quasi-particle degrees of freedom from a high 
density region where they coincide with flavor eigenstates to the vacuum where they coincide with 
mass eigenstates. We write down the most general consistency equations which should be satisfied 
by the effective one-body eigenstates and show that they reduce to the spectral split consistency 
equations for the appropriate initial conditions. 

PACS numbers: 14.60.Pq, 26.30.-k, 02.30.1k 
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1. INTRODUCTION 

Collective oscillations of neutrinos is the process by which neutrino-neutrino scattering contributes to the flavor 
evolution of a sufficiently dense neutrino gas in a way which is somewhat similar to the well known matter-enhanced 
neutrino oscillations in the Sun. Such high neutrino densities are believed to be achieved in the core-collapse supernovae 
in the Early Universe 043; an d possibly in other astrophysical sites. The contribution of the neutrino- neutrino 
scattering to the flavor evolution differs from that of ordinary matter scattering due to the fact that in the former 
the scattered and the scattering particles are of the same kind, giving rise to exchange type forward scattering terms 
[Tol fll| . These terms couple the flavor evolutions of neutrinos with different energies and turn the study of the system 
into a nonlinear many-body problem. A mean-field type effective one-particle approximation to this problem was 
proposed in Refs. [Tof, El an d has been widely adopted in the subsequent studies. However, it became evident that 
analytical solutions of the resulting non-linear evolution equations were needed in order to explore the full range of 
physical scenarios. Recent analytical studies of various special cases, such as the limit in which neutrino-neutrino 
scattering potential becomes dominant, revealed several situations in which neutrinos with different energies oscillate 
collectively 

An algebraic approach to the problem was worked out in Ref. [21| from a many-body point of view. In this approach 
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the flavor evolution of the many-body system is formulated as an SU(2) or SU(3) coherent state path integral for two 
or three flavors, respectively. The evolution operator for the entire system is calculated using both the saddle point 
and the operator product linearization approximations and these two methods were shown to yield the same answer. 
Such an approach is useful in providing a framework to look for exact solutions or systematic approximations. 

There is an increasingly growing literature studying the collective neutrino oscillations (see, e.g., Refs. [22| through 



471). There are, however, several excellent recent reviews that may serve as a starting point in exploring this literature 

HlHi. 



The algebraic approach proposed in Ref. [2l[ is helpful in exploring the hidden symmetries of the system. Hamil- 
tonian describing collective neutrino oscillations possesses an SU(N) / rotation symmetry in the neutrino flavor space 
[2TI [25l [26| . Various collective modes, includin g sp ectral swappings or splittings arise from this symmetry even in 
the inhomogeneous or anisotropic environments |25| . One expects that such a complex nonlinear system may exhibit 
further symmetries. Indeed, several authors noted the presence of various conserved quantities in collective neutrino 
oscillations More recently, it was shown that collective oscillations that maintain coherence can be classified 

by a number of linearly- independent functions [47l ] , implying that scalar products of a unique linear combination of 
the original polarization vectors are conserved. The goal of this paper is to further explore symmetries and conserved 
quantities associated with the collective neutrino oscillation Hamiltonian. 

Our study of the symmetries of the neutrino Hamiltonian is based on the observation that the neutrino-neutrino 
forward scattering Hamiltonian has the form of a spin-exchange interaction as was pointed out by many authors 
earlier. Here, the spin does not refer to the intrinsic spin of the neutrino but to the so called neutrino isospin which 
is defined by introducing a multiplct of neutrino states. Interactions of this type are also encountered in many-body 
systems with pair coupling where the role of the spin is played by the so called quasi-spin. Examples include the 
residual pairing interaction between nucleons in the nuclear shell model and the pairing of valance electrons in the 
BCS theory of superconductivity [48j]. The fact that the pairing Hamiltonian is exactly solvable, which hints at the 
existence of symmetries and associated constants of motion, was pointed out as early as 1963 by Richardso n | 49|| . 
These constants of motion were later identified by Gaudin (50, IHH and others [52| (for a review, see Refs. [H, |54||). 

Here, we consider both the exact many-body Hamiltonian and the effective one-body Hamiltonian describing collec- 
tive oscillations of neutrinos. We formulate the effective one-body picture in terms of the random phase approximation 
(RPA) method whereby the intrinsic consistency requirements are manifested as equations of motion of the neutrino 
polarization vectors. Our analysis includes both the vacuum oscillations and the self interactions of neutrinos. We 
show that, under the single angle approximation, both the exact many-body and the effective one-particle pictures 
possess many constants of motion which were not so far carefully studied. We express these constants of motion in 
terms of the neutrino isospin operators in the case of the exact many-body Hamiltonian and in terms of the neutrino 
polarization vectors in the case of the RPA evolution. In both cases, the constants of motion depend on the param- 
eter which couples the flavor evolution of neutrinos with different energies manifesting the existence of associated 
dynamical symmetries. 

We give analytical expressions for the exact eigenstates and energy eigenvalues of both the many-body Hamiltonian 
and the effective one-particle Hamiltonian. To achieve this goal, we use the method of Bethe ansatz [55] . For the RPA 
Hamiltonian we use the method of Bogoliubov transformations to bring the system into a diagonal form in terms of 
noninteracting quasi-particle states. We show that the quasi-particle picture is useful in offering a formal and intuitive 
description of the spectral splits which were reported to occur in various numerical simulations when the neutrinos 
adiabatically evolve from a region of high neutrino density to the vacuum. Our explanation is complementary to the 
one offered in Ref. 20] in terms of the neutrino polarization vectors. A preliminary account of our results was given 
in 0. 

The organization of this paper is as follow: In Section [3J we briefly review the concept of neutrino isospin and write 
down the many-body Hamiltonian describing the collective neutrino oscillations in terms of the isospin operators. 
We also present the many-body constants of motion in terms of the isospin operators in this section. In Section [J2 
we present the exact many-body eigenstates and eigenvalues which are found in an analytical way using the method 
of Bethe ansatz. In Section 2J we apply the method of RPA to bring the many-body neutrino Hamiltonian into 
its effective one-particle form and we briefly review how the RPA consistency requirements yield the time evolution 
equations of the neutrino polarization vectors. We write down the constants of motion of the RPA evolution in terms 
of the polarization vectors. In this section, we also consider the eigenstates and eigenvalues of the RPA Hamiltonian 
and write down the RPA consistency equations for these eigenstates. At the end of this section we give a brief 
interpretation of the spectral splits in terms of the adiabatic evolution of quasi-particles. The first four sections of 
this paper deal only with neutrinos in order to keep the formulas simple and emphasize the underlying physics. We 
include the antineutrinos in Section [5] and describe how the results of the earlier sections, including the invariants, 
are generalized in this case. We conclude the paper in Section [S] by elaborating the connection between our invariants 
and other invariants described in the literature. 
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2. THE ISOSPIN FORMULATION OF THE PROBLEM 
A. Mass and Flavor Isospin Operators 

In this paper, we consider the mixing between two neutrino flavors. Without loss of generality, we can take one of 
these flavors to be v e and the other to be an orthogonal flavor state that we denote by v x . In other words, v x can be 
either or v T or a normalized combination of them. We denote the fermion operator for an a flavor neutrino with 
momentum p by a a (p) where a = e,x l . The global rotation 

a e (p) — cos0 ai(p) + sin0O2(p) (2.1) 
a x (p) = - sin# ai(p) + cos 9 a 2 (p) 

relates them to the fermion operators dj (p) for the corresponding mass eigenstates for i = 1,2. 
Following the earlier literature, we introduce the flavor isospin operators 

Jp = 4(pK(p) , J P = 4(p)a e (p) , Jp = - (oj(p)ae(p) - 4(p)a x (p)) , (2.2) 
which obey the usual SU(2) commutation relations 

[Jp, J q ] = 2(5 pq Jp , [J p , Jq] — ±SpqJp , (2-3) 

such that we have as many orthogonal SU(2) flavor isospin algebras as the number of neutrinos. It follows from the 
definitions given in Eq. (|2.2j) that each flavor isospin algebra is realized in the spin-1/2 representation and that the 
electron neutrino is taken to be isospin up. We note, however, that sometimes the opposite convention for the isospin 
doublet is used (see, for example, Ref. [21j]V 

We introduce the following summation convention for all quantities labeled by neutrino momentum p: 

a u = J2 a p and A = J2 A ^ • ( 2 - 4 ) 

|p|=p w 
Here ui is the vacuum oscillation frequency for a neutrino with energy p. It is given by 

U - f (,5> 

with 5m 2 = m\ — m\ where m, is the mass of the state Vi. In the case of flavor isospin operators, for example, J u 
represents the total flavor isospin operator of all neutrinos with the same vacuum oscillation frequency uj whereas 
J represents the total flavor isospin operator of all neutrinos. Since the operators J u and J are sums of individual 
SU(2) operators, their components also obey the SU(2) commutation relations. The corresponding total flavor isospin 
quantum numbers can take several values. Note that in Eq. (|2.4I) and in all subsequent equations sums are to be 
taken over all occupied neutrino states. Also note that we use boldface letters to indicate vectors in momentum space 
(e.g. p) and arrows to indicate vectors in flavor space (e.g. J). Throughout this paper, we use both the Cartesian 
basis (x, y, z) and the cylindrical basis (+, — , z) for vectors in flavor space. These two bases are related by 

J±=J X ± iJy . (2.6) 

One can similarly introduce the mass isospin operators 

J+ = oJ(p)a a (p) , Jp - 4(pK(p) > J p=1 ( at i(p)ai(p) - 4(pMp)) , (2.7) 

which also obey the SU (2) commutation relations given in Eq. (|2.3[) . As can be seen from the definitions given in 
Eq. (|2.7j) . we take the mass eigenstate vi to be isospin up. The total mass isospin operators are defined as in Eq. 
(El. 



1 In general, additional quantum numbers are needed in order to distinguish neutrinos with the same momentum. But to keep our 
formulas simple we do not explicitly include additional quantum numbers in our notation. Instead, one can view p as a multiple index 
like (p, si, S2, ■ ■ ■ )• 
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Eqs. (|2.1I) imply that the particle operators in the mass and flavor bases are related by a unitary transformation. 
The most general unitary transformation operator in flavor space can be written as 

Up = e E P *vJ+ e E P in{i+\z P \ 2 )j; e - E P z P J P (2.8) 
where z p — e tSp tan# p is a complex parameter. The effect of such a general transformation on particle operators is 

J7p t Oi(p)f7 p = cos6» p ai (p) - e iS <> sin9 p a 2 (p) (2.9) 
U p ^a 2 (p)U p = e~ t5 v sin6> p oi(p) + cos6» p a 2 (p) 

i.e., each momentum mode is rotated by a different angle and acquires a different (relative) phase. The transformation 
from mass to flavor basis given in Eqs. (12.1[) is a special case of Eqs. (|2.9[) for which the rotation angle is equal to the 
vacuum mixing angle and the phase is zero for all momentum modes. In other words, we can write 

o e (p) = UUx{v)U and a x (p) = U^a 2 (p)U , (2.10) 

where 

U = e zJ+ e Ki+W a )^ e - zJ ~ (2.11) 

with z p = z = tan 9 for all p. Note that we use the subscript p in denoting the operator U p defined in Eq. (|2.8I) in 
order to emphasize that it imposes a different rotation on each momentum mode. Although the operator U is only a 
special case of U pi it carries no such index because it induces the same transformation on all momentum modes. Also 
note that the rotation parameter z p is not subject the summation rule introduced in Eq. (|2.4j) . Instead, its index 
simply indicates its dependence on momentum. 

The operator U p given in Eq. (|2.8[) also induces a transformation on the isospin operators as follows: 

UlJpU P = ™s29 p j; + ie**- sm20 p J+ + sm26 p J~ , 

UlJfU p = cos 2 9 P J+ - e~ iS » sin 29 P J£ - sin 2 9 p J p , (2.12) 

U 1^p U p = cos2 pJp - ^ s ' m 20pJ P - e 2lS * sin 2 6 p J+ . 
In particular, mass and flavor isospin operators are related by 

J P = U^J P U (2.13) 



which implies 



J* = cos 29 Jp 2 + sin 26> + Jp 
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J+ = cos 2 9J+ - sin 29 J* - sin 2 9J p , (2.14) 
J p = cos 2 6J- - sin 29 - sin 2 9J+ . 

The inverse transformation from flavor to mass basis can be found by substituting —9 in place of 9 in Eq. (|2.14p . Note 
that Eqs. (|2.11[) and (j2 . 13[) follow from the fact that the total isospin operator J is the generator of global rotations 
in flavor space. 



B. The Hamiltonian and the Quantum Invariants 

The many-body Hamiltonian describing the vacuum oscillations of a group of neutrinos is given by 

,2 m 2 



H v = J2[ ^i(p)«i(p) + 2^4(P)«2(P) ) ' ( 2 - 15 ) 



One can write this Hamiltonian in terms of the neutrino isospin operators defined in Section 12 AI First note that 
since the neutrinos only exchange their momenta with the forward scattering, the total number of neutrinos in each 
momentum mode is constant. As a result, the term 

E ! ^4^(a t 1 (p)«i( P ) + 4(P)« 2 (P)) (2.16) 
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is proportional to identity for a given number of particles and can be subtracted from the Hamiltonian without any 
consequences. Subtracting this term from Eq. (|2.15[) and using the definitions given in Eqs. (|2.5p and (|2.7[) together 
with the transformation given in Eqs. (|2.14|) one finds 



Sm 2 
~2p 



Here B is the unit vector which points in the mass direction. Its components are given by 

B = (0, 0, -l) m „. = (sin 29, 0, - cos 20) flavor (2.18) 

in mass and flavor bases, respectively. 

In writing the Hamiltonian which describes the self refraction of neutrinos, one should take into account the fact that 
there are two kinds of neutrino-neutrino scattering diagrams which add up coherently during the neutrino propagation. 
One of them is the forward scattering diagram in which there is no momentum transfer and both neutrinos remain in 
their original states after the scattering. These diagrams give rise to a diagonal refraction potential in the flavor basis, 
similar to the MSW potential. The other one is the exchange diagram in which the neutrinos exchange their states 
after the scattering. These diagrams give rise to a non-diagonal refraction potential in the flavor basis which is the 
source of the nonlinearity of the neutrino self refraction problem [13, EH Ej] . Taking into account the contribution 
of both kinds of diagrams mentioned above, the effective many-body Hamiltonian describing the self refraction of a 
dense neutrino gas can be written as [56j | 

H vu = — =^^^(1 -cosi? pq )[ 4(p)a e (p)4(q)a e (q) +at(p)a a: (p)4(q)a ;E (q) 

p q 

+4(pK(p)4(qK(q) + 4(pK(p)<4(qK(q)) ] • (2.19) 

Here V is the quantization volume and $ pq is the angle between the momentum modes p and q. The factor 1 — cos i9 pq 
guarantees that those neutrinos traveling in the same direction do not undergo scattering. 

One can write the neutrino self refraction Hamiltonian in terms of the neutrino isospin operators defined in Section 
12 Al as follows: 

ff w = ^^£(l-cosi? pq ).V V (2-20) 

Here, Eq. (12.20[) can be directly obtained from Eq. (12.19[) by using the definitions given in Eq. (|2.2p and discarding 
those terms proportional to identity. Eq. (|2.20j) tells us that in the language of flavor isospins the neutrino-neutrino 
interaction takes the form of a spin exchange interaction. This is an expected result because the effective Hamiltonian 
given in Eq. (|2.19[) consists only of those neutrino-neutrino interactions in which neutrinos either keep or exchange 
their momenta. 

The flavor evolution of a dense neutrino gas (in the absence of any other background) is described by the sum of 
the vacuum oscillation term and self interaction term 



ff = ^o;£- J w +/x^(l-cos<!? pq ) J p - J q . (2.21) 

UJ 

Here we defined 



p.q 



y2Gw , . 

,i = . (2.22) 

In this study, we adopt what is commonly referred to as the "single angle approximation," i.e., we will assume that the 
term involving cos $ pq in the Hamiltonian averages to zero so that the neutrinos traveling in different directions which 
are otherwise identical undergo the same flavor evolution. With this simplification, the total Hamiltonian becomes 

H = -^ujX + iiJ-J (2.23a) 



Y^uB ■ J u + fj,J- J . (2.23b) 
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Note that the self interaction term has the same form in both mass and flavor bases because these two bases are 
related by a global rotation as described in Eq. (I2.13[) which leaves all scaler products invariant. 

It has been pointed out by many authors that the Hamiltonian given in Eq. (|2.23j) is analogous to the Hamiltonian 
of an interacting spin system, i.e., a group of spins interacting with a position dependent external magnetic field and 
with each other via spin exchange interaction. In the mass basis the external magnetic field points in the —z direction 
whereas in flavor basis it points in the direction of the unit vector B = (sin 26, 0, — cos 26). This analogy makes it 
easier to see that the length of each isospin is conserved: 

L U = J U -J U [H,L u }=0. (2.24) 

Similarly, the total isospin component in the direction of the external magnetic field is also conserved: 

Co = B-J [H, C ] = . (2.25) 

It is worth noting that the Hamiltonian in Eq. (I2.23j) also appears in connection with the pairing problem. For 
example, in the BCS theory of superconductivity [48J the (reduced) pairing interaction between the valance electrons 
is described by the Hamiltonian 2 

H BCS = J2 - GT+T- . (2.26) 

k 

Here G is a constant which represents the strength of the pairing interaction. Within the context of the BCS model, 
it is assumed that G > leading to an attractive pairing and making the formation of Cooper pairs energetically 
favorable, e k are the degenerate single particle energy levels which can be occupied by pairs of spin-up and spin-down 
electrons (i.e., the Cooper pairs). The operators 

** = c lr c U i *k = c ^Cfet and t% = i (c^Cfet + c^c^ - l) (2.27) 

are called quasi-spin operators and they also obey the same SU(2) commutation relations as given in Eq. (|2.3[) . In 
the quasi-spin scheme, a single particle state k has quasi-spin up if it is occupied by a pair and quasi-spin down if it 
is not. The operator 

T = J2*k (2-28) 

k 

represents the total quasi-spin. As can be easily verified, the third component of the total quasi-spin is a constant of 
motion, i.e., 

[H BCS ,T*]=Q, (2.29) 

which is analogous to the conservation of Co mentioned in Eq. (|2.25[) . From Eq. (|2.27p we see that 2T Z is equal to 
the total number of electron pairs in the system. 

The Hamiltonians given in Eqs. ()2.23a[) and (|2.26l) are very similar. In fact, if we make the substitutions uj — > 2e k 
and J7p — > tk and keep in mind that the 'missing' term GT Z (T Z — 1) in Eq. (|2.26[) is proportional to identity for a given 
number of pairs (i.e., has no effect on the evolution), then we see that the two models are mathematically identical 
up to an overall minus sign (remember that \x is by definition a positive quantity) . The correspondence between these 
two models can be summarized as follows: Single particle states with energy e k in the BCS model correspond to the 
neutrino states with oscillation frequency uj. Both sets of states are multiply degenerate. A state ek which is occupied 
(respectively, unoccupied) by a pair in the BCS model corresponds to a neutrino state with oscillation frequency uj 
occupied by a neutrino in V\ (respectively, v<i) mass eigenstate in the neutrino model. 



2 The most general (i.e., non reduced) pairing Hamiltonian has the form 

-Hbcs = ^2 2tktl — G ^2 c fefe'*fc * fc / i 
fc k,k' 

where the dimcnsionlcss coefficients c k y lead to a state dependent pairing strength similar to the coefficient 1 — cos i? pq in the neutrino 
Hamiltonian. 
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It was first shown by Richardson in 1963 that the BCS Hamiltonian given in Eq. ( 2.261) c an be diagonalized 
analytically with the method of algebraic Bethe ansatz [4!| [55[ . Later work by Gaudin 50l l5l| and others [52l - [54j 
revealed that the integrability of the Hamiltonian in Eq. (|2.26j) derives from the existence of a set of constants of 
motion or quantum invariants. These results can be easily carried over from BCS model to the neutrino model. For 
example, it can be shown that the neutrino Hamiltonian given in Eq. (|2 . 23[) has the following constants of motion: 

h u = B-J u + 2fi V J -^L . (2.30) 

The operators given in Eq. (|2.30p are known as the Gaudin magnet Hamiltonians. It is straightforward to show that 
these operators commute with one another and with the Hamiltonian given in Eq. (I2.23p . i.e., 

[h u ,h u >}=0 and [H,h u ]=Q (2.31) 

is satisfied for every u> and co' . Note that the invariants given in Eq. (12.30[) are independent from one another and 
from the invariants mentioned in Eq. (|2.24p . However, the invariant Co mentioned in Eq. (|2.25p is the sum of the 
Gaudin magnet Hamiltonians 

C = J2 h " ( 2 - 32 ) 

Ul 

and thus is not an independent invariant. The Hamiltonian itself is also a linear combination of the invariants given 
in Eqs. (g^D and 

Ul Ul 

This tells us that if the neutrinos occupy Q different energy modes, then the system has 2il independent constants of 
motion given in Eqs. (|2.24p and (|2.30p . These invariants can be expressed in several alternative forms which may be 
useful in different applications. In Appendix \^ we present two alternative ways of writing down the invariants. 
From the definitions given in Eq. (|2.7p . one can see that the constant Co can be written as 

Co = -J z = 2 2 1 (2.34) 

where 

& =H a i(p) a i(p) and N 2 = ^24(p)a 2 {p) (2.35) 

p p 

denote the total particle number operators for the mass eigenstates v\ and 1/2, respectively. Note that in this paper 
we denote the particle number operators with hats (e.g., Ni) and the corresponding particle numbers without hats 
(e.g., Ni). Together with the fact that the total number of neutrinos N = Ni + N2 is constant, the conservation of 
Co mentioned in Eq. (|2.25p ensures that the Hamiltonian preserves the number of neutrinos in each mass eigcnstatc, 
i.e., we have 



H,Ni 



, (2.36) 



for i = l,2. 

Although Co, which is a combination of the invariants h u , can be simply expressed in terms of particle number 
operators as in Eq. (I2.34p . this is an exceptional situation. Apart from this particular case, the invariants h u or their 
combinations cannot be written in terms of the particle number operators. This is evident from the fact that the 
terms J w • in Eq. (|2.30l) are diagonal in neither the mass basis nor the flavor basis and that these terms disappear 
only in the particular combination given in Eq. (12 .32[) . 



3. DIAGONALIZATION OF THE HAMILTONIAN 



The existence of constants of motion is a manifestation of exact solvability of the pairing Hamiltonian given in Eq. 
(|2.23p . As mentioned earlier in previous sections, the pioneering work in this direction was that of Richardson who 
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showed that the application of the Bethe ansatz method to the pairing Hamiltonian yields its exact eigenstates and 
eigenvalues in an analytical way. This exact solvability of the pairing Hamiltonian has been studied and exploited 
extensively thereafter both in the context of the nuclear shell model and in the context of BCS theory. In this section, 
we will review this procedure for a self interacting neutrino gas. 

The Bethe ansatz method, if applicable to a problem, usually reveals a general and simple functional form for 
the many-body eigenstates. In the case of the self interacting neutrino gas, the many-body eigenstates turn out to 
resemble those of an harmonic oscillator in that they can be obtained by repeated application of a step operator on a 
lowest weight state. The crucial difference is that in the neutrino case the step operator is parametrized by a complex 
number. This parameter takes on different values at each step and should be determined by solving the equations 
of Bethe ansatz as described in more detail below. The number of equations that one is required to solve in order 
to find the full spectrum exactly is of the order of the number of particles. For a dense neutrino gas, this is clearly 
unfeasible. However, if an approximate method of determining the Bethe ansatz variables can be found, then the 
simple step operator form of the eigenstates allows one to write down an (approximate) evolution operator for the 
system. Alternatively, one can consider a relatively small number of neutrinos occupying a very small volume so as to 
yield a large density. Mimicking the conditions in a dense environment in this way allows one to explore many-body 
physics with the exact eigenstates. One can then increase the number of neutrinos and the volume that they occupy 
keeping the density constant. The resulting limit of the Bethe ansatz equations have been studied extensively in the 
context of the BCS model (see, for example, Refs. [H, [57l - [6l| and the references therein) and in connection with 
matrix models that appear in 2D gravity [62]. This limit will be reviewed in Section [3 Bl 

We would like to note that a study of the general functional form of the eigenstates as revealed by the Bethe ansatz 
method can itself give insight into the collective behavior of the system. For example, it was shown in Ref. [g3[ 
that the expectation value of the step operator in the effective one-particle approximation to the pairing problem 
is the generating function for the canonical variables of the system and that the exact analytical solutions of the 
corresponding RPA equations of motion can be obtained in terms of these variables. 

A. The Bethe Ansatz Method 

The eigenstates of the neutrino Hamiltonian given in Eq. (I2.23P can be easily found in the two opposite limits of 
the parameter fi defined in Eq. (|2.22[) . 

• As ji approaches to zero, neutrinos occupy a larger and larger volume and neutrino-neutrino scattering 
becomes negligible. In this limit, the Hamiltonian consists only of vacuum oscillations given in Eq. (|2.17l) . The 
eigenstates in this limit are simply the tensor products of mass eigenstates for individual neutrinos. Speaking 
in terms of the interacting spin system analogy, in this limit the spins interact only with the external magnetic 
field which is in the — z direction in mass basis. The eigenstates are those in which each spin is either aligned 
or anti-aligned with the external magnetic field, i.e., 

\U\ V\ V\ . . .), \U 2 V\ V\ . . .), \l>i V2 V\ ...), \V\V\V-z ...), ■■■ (3.1) 

• As /i approaches to infinity, neutrinos are crowded into a smaller and smaller volume and eventually the 
neutrino-neutrino scattering term becomes dominant. In this limit, one can ignore the vacuum oscillations and 
write the Hamiltonian as 

ffoo = lim H = fif-J. (3.2) 

The total isospin quantum number j is a scalar quantity so that we can use it without any references to mass or 
flavor bases. Note that j can take several values starting from or 1/2 (depending on whether we have an even 
or odd number of neutrinos) up to j max = N/2 where N is the total number of neutrinos. In the limit where 
fi — > oo, both the total mass isospin states \j, m) m and the total flavor isospin states \j, m)/ are eigenstates of 
the Hamiltonian with the same energy, i.e., 

Hoo\j, m) m = + m) m and Hoo|j, m) f = + m) f . (3.3) 

It is clear from Eq. (j2.13[) that these two sets of states are related by 

\j,m) f = U^\j,m) m . (3.4) 

Since the operator U involves only the total isospin operators (see Eq. (12. 1 1[) ^ it cannot change the representation 
in which a state lives, i.e., it preserves the value of j in accordance with the above arguments. 
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To illustrate the use of Eq. (|3.4[) . let us consider a state in which all neutrinos are v e . Since this corresponds to 
having all flavor isospins up, this state is the highest weight state of the total flavor isospin algebra, i.e., 



■ ) = \u^uJ)f- (3.5) 
This state can be converted to mass basis using Eqs. (|2.1ip and (|3 .4[) . The result is 

\3*~,3*J)f = tf\j m ,j^)m (3.6) 

= 5 (cosdy™ - m ( S mOy~» +m J T . (2j "" )! |j m „,m) m . 

m= T/ max V Wm« - m) ! (j max + to) ! 

We see that the state |j max ,j max ) / is a linear combination of the states |j max ,TO,) m , which all live in the j = j max 
representation and have the same energy in the limit where fx — > oo. 

We next present the eigenstates and eigenvalues of the collective neutrino oscillation Hamiltonian given in Eq. 
(|2.23[) away from those two limits, i.e., for < \i < oo. We begin by noting that Eq. (|2.36|) implies that all eigenstates 
of the collective oscillation Hamiltonian are also eigenstates of the number operators N\ and iVa, i.e., each many-body 
eigenstate has a definite number of neutrinos in v\ and v% states. 

In general, the total isospin states \j,m) m and \j,m)f mentioned in Eq. (|3.3p are no longer eigenstates of the 
Hamiltonian for a finite value of fi but the highest and lowest weight states \j, ±j) m of the total mass isospin continue 
to be eigenstates. Let us consider, for example, the states in which all neutrinos occupy the same mass eigenstate, 
i.e., they are all v\ or all v 2 . Speaking in the language of isospin, this is equivalent to having all mass isospins up or 
all mass isospins down, respectively. These situations respectively correspond to the total isospin states |j max , j max ) m 
and |j max , -j max ) m : 

\u^j^)m = l[a\(p) |0) and | Jmax , - Jmax ) m = Y[4(p) |0) . (3.7) 
p p 

One can easily show that these states are eigenstates of the Hamiltonian given in Eq. (|2.23p with the respective 
energies 

UJ UJ 

Here we denote the total number of neutrinos in an individual energy mode by N u . The other highest and lowest weight 
states \j, ±j) m of the total mass isospin are also eigenstates of the Hamiltonian. To see this, let us first consider the 
total isospin quantum number j u of an energy mode which can take several values from or 1/2 (depending whether 
we have an even or odd number of neutrinos in this mode) up to N u /2, possibly with many multiplicities. For each 
energy mode, we can write 



±juj\juj,±juj) m ■ (3.9) 



Therefore the states 



(3.10) 

lii ~j)m = b'wii ~jux)m ® \ju21 ~joj 2 )m ® ' ' ' \jwni ~3osn)m 

are simultaneous eigenstates of all 3* (see Eq. (|3. 91) ) as well as the quadratic operator ^7 • (see Eq. p.3p ). As a 
result, the states given in Eqs. (|3.10p are eigenstates of the Hamiltonian (|2.23[) for all values of fi with the respective 
energies 

E {+3) uju + M {j + 1) and E { _ 3) = ujj^ + fj,j (j + 1) . (3.11) 

In Eqs. (|3.10p , we used to denote the number of different energy modes. In this equation the total angular 
momentum j = j Pl + j P2 + • • • +j PU may come with many multiplicities because in general more than one combination 
of j Pl ,j P2 , ■ ■ ■ ,jpn may correspond to the same sum j. Such multiplicities are inherent in the addition of angular 
momenta and one usually introduces additional quantum numbers to distinguish the resulting degenerate total angular 
momentum states. However, in this paper we avoid introducing such quantum numbers to keep the formulas readable. 
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Since the energies given in Eq. (|3.11[) depend only on j u and j, this omission does not lead to confusion in what 
follows. Finally note that Eqs. (|3.7[) and (I3.8P correspond the unique special case of Eqs. p.lOj) and p. lip in which 
we have j u = N^/2 for all energy modes and thus j = j max . 

The eigenstates of the Hamiltonian other than the highest and lowest weight states \j, ±j) m can be found with 
the Bethe ansatz technique. This method is based on writing down a trial state depending on a set of unknown 
parameters which are called Bethe ansatz variables. The values of these parameters are determined subject to the 
requirement that the state be an eigenstate of the Hamiltonian. For the particular problem at hand, the trial state is 
formed with the help of the Gaudin operators 

s ± (0 = E^7- ( 3 - 12 ) 

LI — f 

Here £ is a complex parameter which will later play the role of a Bethe ansatz variable. In order to demonstrate the 
method of Bethe ansatz, let us consider the state \j, —j) m defined in Eq. (|3.10|) . This state contains Ni = N/2 — j 
neutrinos in mass eigenstate v\ and Ni = N/2 + j neutrinos in mass eigenstate v-i- The operator Q + (£) turns one v% 
into V\ such that the state 

\0 = Q + (0\j,-j)m (3.13) 

contains one more v\ and one less v-i. It is, in fact, a linear superposition of many such states depending on the 
unknown parameter £. In order to find for which value(s) of £ the state in Eq. (|3.13l) is an eigenstate, we act on it 
with the Hamiltonian given in Eq. (|2.23l) . The result is 

HQ + (0\j, -j) m = - £ - 2/y) C+(0|j, -J) m + Z^l) J+lj > ~ j)m ■ (3 - 14) 

Here E(-j) is the energy of the state \j,—j) m given in Eq. (13. lip . In order for the state Q + (£)|i, — j) m to be an 
eigenstate of the Hamiltonian, we should choose £ in such a way that the second term on the right hand side of Eq. 
p,14p vanishes. This tells us that £ should obey 



E 



- 1 (3.15) 



u) — £ 2/i 



Eq. (|3.15p is called a Bethe ansatz equation and in general it has several solutions. It is clear from Eq. (|3.14p that 
every £ which satisfies Eq. (|3.15p gives us an eigenstate in the form of Eq. (|3.13p with the energy 

E(0 = E { _fi - i - 2nj . (3.16) 

One can extend this line of thought to find more eigenstates of the Hamiltonian. In the most general case, the 
Gaudin operator plays the role of a one parameter step operator, i.e., the eigenstates can be obtained by its repeated 
application on a lowest weight state but the parameter £ takes on different values at each step. It can be shown that 
a state in the form 

&,...£«) = Q + (£i)Q + (6)--.Q + (£«)li,-i) m (3.17) 

is an eigenstate of the Hamiltonian with the energy 

Efa, Z K )=E i _ j) -J2Z*-w(2j-K+l) , (3.18) 

a=l 

if the Bethe ansatz variables £i,£2j • • • obey the Bethe ansatz equations 

Eqs. (13.191) form a set of k complex equations in n complex variables that have to be simultaneously satisfied for 
every a = 1, 2, . . . , k. As stated earlier, the state \j, —j) m contains respectively Ni = N/2 — j and N 2 = N/2 + j 
neutrinos in mass eigenstates U\ and i>2- Since each one of the Gaudin operators Q + (£ Q ) transform one v 2 into U\, 
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corresponding occupancies for the state in Eq. (|3.17l) are Ni = N/2 — j + k and N 2 = N/2+j — k. It should be noted 
that if we add one more Gaudin operator to the state in Eq. (|3.17[) to transform one more neutrino from v 2 to Ui, 
then the Bethe ansatz equations ()3.19j) are modified and become a new set of k + 1 equations in k + 1 variables. This 
new set of equations are different from the earlier ones with k variables and have different solutions which could be 
denoted as (£1, £2, £k+i)- Each time one adds one more Gaudin operator to Eq. Q3.17P to find new eigenstates 

with different occupancies, one has to solve a new set of coupled algebraic equations. Note that the example described 
in Eqs. (13.13113.1^1) is a special case of the general scheme described in Eqs. (|3.17M3.1§)) with k = 1. In this particular 
case, the sum on the right hand side of Eq. (|3.19p vanishes because there is only one Bethe ansatz variable. 

It is obvious from the above remarks that it becomes more and more difficult to find eigenstates as we flip more and 
more neutrinos from v 2 into V\. However, a symmetry transformation between the eigenstates of the Hamiltonian is 
helpful in reducing the number of Bethe ansatz equations that has to be solved. In order to present this symmetry, 
let us first observe that the operator 

T = e~™ J * (3.20) 

transforms v\ and v 2 neutrinos into each other. In particular, it exchanges the states \j, —j) m and \j,j) m defined in 
Eqs. (|3~T0f 

T%-j)m = \3,3)m (3.21) 

and it also transforms the isospin operators as 

T*J±T = J2 and T^J*T = — J* . (3.22) 
As a result, the Hamiltonian given in Eq. (|2.23ap is transformed as 

T*HT = H' = ^uJZ + J '■ J . (3.23) 

Elementary quantum mechanics tells us that if a state \ip') is an eigenstate of the Hamiltonian H' — T'HT with a 
particular energy, then the state \tp) — T\ip') is an eigenstate of the Hamiltonian H with the same energy. Therefore, 
if we find the eigenstates and eigenvalues of the Hamiltonian H' with the method presented above (by making the 
change uj — ¥ — u) and then transform the results with the operator T, we arrive at the eigenstates of the Hamiltonian 
H given in (|2.23[) . In this way, one can show that the states 

\&,&,...a = Q-&)Q-(Q...Q-(0\j,j) m (3.24) 

are eigenstates of the Hamiltonian H with the energy 

E(^,^,...,Q=E (+j) +'£C-K t ,(2j- K + l) , (3.25) 

a=l 

if the variables £[, £ 2 , . . . , obey the Bethe ansatz equations 



These Bethe ansatz equations are similar to the ones given in Eq. (|3.19j) except that the sign of the constant term 
is reversed. The term which appears in Eq. f|3 . 25[) is the energy of the state \j,j) m given in Eq. (|3.11l) . The 

state \j,j) m which is defined in Eq. (13.10[) contains respectively Ni = N/2 + j and N2 — N/2 — j neutrinos in mass 
eigenstates v\ and Each Gaudin operator Q~{£/ a ) transforms one v\ into v 2 so that the corresponding occupancies 
of the state in Eq. (|3.24[) are N\ = N/2 + j — n and N 2 = N/2 — j + n. We see that the occupancies of the mass 
eigenstates v\ and v 2 are reversed in the state given in Eq. (|3. 241) as compared to the state given in Eq. (|3.17|) 
but both states are reached by solving the same number of equations. Therefore, depending on the occupancies of 
the mass eigenstates v\ and v 2 , either the method presented in Eqs. (|3.1 7H3.19")) or the method presented in Eqs. 
is more economical to use. 



12 



+1 



E=- 



, , +lo c 

S7 



~J«, ~./«,, ~J<« 
• • — : — • — 



2fi 



• — *x 



C0 2 W 3 



-lo r 



+1(1 



?S +1( 



5 8 +U 



FIG. 1: The Bethe ansatz equations (|3.19|) can be viewed as the stability conditions for a group of point charges in two 
dimensions as shown here. 



B. Electrostatic Analogy 

The Bethe ansatz equations given in Eqs. Q3.19P and (|3 . 26[) can be viewed as the stability conditions for a group 
of point charges in two dimensions 3 [59(. Let us consider Eqs. ()3.19p for which the electrostatic model is depicted in 
Fig. Q] There are two kinds of charges in this scheme: the fixed negative charges and the free positive charges. For 
each energy mode we have a fixed point charge of magnitude — j u pinned at the point uj on the real axis and for each 
Bethe ansatz variable we have a free charge of magnitude +1 whose equilibrium position yields £ Q . There is also a 
constant and uniform electric field in the —x direction as shown in the figure. 

In two dimensions, one can use a complex coordinate z — x + iy instead of two real coordinates (x,y). The 
electrostatic potential at the point z created by a charge q sitting at the point zq is proportional to qbi \ z — zq\. As a 
result, the total electrostatic energy of the charge configuration shown in Fig. [T| is proportional to 

V(xj-J2R<^)-j-J2^ Re ^-l E J<-iwHw-w'|+I><-H&-w| ■ (3-27) 

a " u) a ,P cj.cj' a.u) 

The free charges come to an equilibrium when the electrostatic potential energy reaches a local minimum. The 
minimum energy condition is obtained by setting 

dV 

a5 = ' <3 ' 28) 

for every a = 1, 2 . . . , k. It is a straightforward calculation to show that Eqs. (|3.27j) and (|3.28j) lead to the Bethe 
ansatz equations (|3.19|) which now tell that the total force on each free charge is zero. A similar electrostatic picture 
for the Bethe ansatz equations (|3.26p can be obtained by reversing the direction of the electric field in Fig. [TJ Note 
that since the solutions of Bethe ansatz equations come in complex conjugate pairs, the organization of the free 
charges is always symmetric with respect to x axis. 

The electrostatic analogy described here provides an intuitive picture and can be particularly useful in finding the 
solutions of Bethe ansatz equations in the limit of a large number of neutrinos. In this scheme, the total amount of 
fixed charge on the x axis is equal to —j. Although j takes all values between and N/2, for most eigenstates its 
value is very large. The total amount of free charges is given by 

« = j + (JVi - N 2 )/2 . (3.29) 
Here we assume that j and K are of same order of magnitude which is true for most eigenstates. 



Alternatively, one can imagine infinite lines of charge in three dimensions perpendicular to a plane. 
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Suppose that we start with a small number of neutrinos in a very small volume so that the neutrino density is in 
the regime where the neutrino-neutrino interactions are important. In this case we have a small amount of charge 
on the x axis, a few free charges for a typical eigenstate and a small external electric field in the — x direction. Note 
that, according to the definition of /i given in Eq. (|2.22p . the external electric field — 1/2/j is proportional to the 
volume occupied by the neutrinos. Now, suppose that we increase the number of neutrinos and the volume which 
they occupy by keeping the density constant. In this case the external electric field grows proportionally while the 
total fixed charge on the x axis and the number of free charges also increase. For a realistic neutrino spectrum, fixed 
charges form a continuous charge distribution on the x axis. However, for the sake of this discussion, we assume that 
they are combined into discreet energy bins such that each j u is of the order of j. Numerical solutions of Bethe ansatz 
equations suggest that in this limit the free charges form (piece-wise) continuous distributions. Such solutions have 
been worked out in the context of the electron pairing in superconductors (53l l57T - l6i| . An intuitive way of dealing 
with this limit was introduced in Ref. [58[ and is based on an expansion of the total electrostatic field in powers of 
1/k. In what follows we closely follow this reference. 

Let us begin by considering the total electrostatic field of the system in Fig. [1] It is given by 

F ^ = -f+E-^ + E -\ ■ ( 3 - 3 °) 

Using the Bethe ansatz equations (|3.19[> . it is straightforward to show that the electrostatic field obeys the following 
differential equation: 

Here H(uj) is the electric field produced by the free charges at the position of the fixed charge — j u . It is given by 

H( U ) = J2 77V ■ ( 3 - 32 ) 

a=l W ~ «° 

H(uS) can also be expressed as a contour integral 

H(u) = ^- I • (3.33) 

2tt* Jc u — z 

Here C is a contour which encloses only those singularities of the field F(z) due to the free charges. The equivalence 
of Eqs. p.32p and (|3.33l) can be shown by direct substitution of the field F(z) given in Eq. (|3.30j) in Eq. (I3.33|) . 
Let us now consider the multipole expansion of the electrostatic field F(z): 

oo 

F(z) = F^z- m . (3.34) 

m— 

A direct expansion of Eq. (]3.30[) gives these multipole moments as 

FV\z) = ~ , (3.35) 
F^{z) = K-j, (3.36) 

K 

F (n) (z) = " n j« n>2, (3.37) 

a — l lo 

for \z\ > max(w, \Ca\)- The zeroth order moment gives the value of the field at infinity and therefore is equal to the 
value of the constant field. The first order moment is equal to the total amount of charge, as expected. By comparing 
Eqs. (|3.37[) and (|3.18j) one can see that the energy of an eigenstate can be found from the second order moment of 
the corresponding electrostatic field. 

Richardson's method is based on solving the differential equation given in Eq. (13.31)) for the leading terms in an 
expansion of the form 



F{z)=Y j F r {z) . (3.38) 
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Here the term F r (z) is assumed to be of the order of k 1 ^ . In the present case, we will consider only the r = 
term which is of the order of k. It is reasonable to assume that those terms on the right hand side of the differential 
equation which are of the form 53 w iw(- ■ ■ ) are of the order of j. Since the volume is increased in proportion with 
the number of neutrinos, we also treat the external field — 1/2 fi as being of the order of j. Finally, Eq. (|3.32[) tells us 
that H u is of the order of k. Taking these into account and substituting the expansion given in Eq. (|3.38[) into the 
differential equation ()3.3ip we find 



where H (uj) is given by 



tf..U., = -L/^. (3.40) 
2ni J c lu — z 



Eq. (|3.39p is an integral equation which includes only the fixed charges as parameters and should be solved to 
find the leading order electrostatic field created by the whole configuration including the free charges in their stable 
configuration. For this reason, it can be thought as a replacement for Bethe ansatz equations. Once a solution of Eq. 
(|3.39j) is obtained, one can determine the leading order contribution to the energy of the corresponding eigenstate 
from the second order moment of the field. The locations of the free charges can be found from the singularities (or 
the branch cuts) of the field other than the fixed charges and these can be used to write down the corresponding 
eigenstate itself. 

However, instead of directly solving Eq. (|3.39|) . a more practical way is to form its solutions on physical grounds. 
The typical way to proceed is to start with a small number of neutrinos in a small volume and see the way the 
solutions of Bethe ansatz equations organize themselves as the number of neutrinos are increased while the density 
is kept at a constant value. As mentioned before, numerical simulations suggest that the free charges coalesce to 
form (piece-wise) continuous charge distributions. One can guess the limiting shape that they will assume possibly in 
terms of some unknown parameters. These parameters can later be determined either by forming a field Fq(z) which 
describes such a distribution to the leading order and substituting it in Eq. Q3.39P or, in simple situations, directly 
from self consistency requirements as described in the following example. 

For example, it is known that a solution of Bethe ansatz equations exists in which all free charges organize themselves 
into a single arc extending from a point a to a point o* (the parameters to be determined). The situation is similar 
to the one depicted in Fig. [T] except that the free charges now form a continuous curve. In the case of the BCS 
model, this solution is known to lead to the BCS ground state. Since the singularities due to the free charges form a 
continuum, the field Fo(z) is expected to have a branch cut along the arc such that 

— / F (z)dz = k (3.41) 
2m J c 

is satisfied for any closed path C enclosing the free charges. Note that it is reasonable to also demand that the field 
F (z) describes the fixed charges and the external field correctly, i.e., 

lim(z - uj)F (z) = j u and lim F (z) = — - (3.42) 

2— >CJ Z— >00 2fl 

because both j u and the external field are of the order of j. In this case, the zeroth and the first order moments given 
in Eqs. (|3.35p and (|3.36p will receive no corrections from the higher order terms in the expansion of electrostatic field 
given in Eq. ([335)) . 

In Ref . [sgj , it is argued that a candidate for the field F (z) with these properties is given by 

F (z) = -[(z-a)(z-a*)] 1/2 J2) 77 \ • ( 3 - 43 ) 



This field is of the order of j and has a branch cut along an arc which extends from a to a* . It also has singularities 
at the points w and a constant limit as z — > oo as required by the above arguments. The first three moments of the 
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field Fo(z) are given by 

U) 1 1 

= £ ¥^ ■ (3-45) 
\u) — a\ 

UJ 1 1 

F o -2. > ( 3 - 46 ) 

where a = A + iA. One can try to determine the value of the parameter a by demanding that Fq(z) reproduces the 
zeroth and the first order moments without any contributions from higher order fields as mentioned above. Equating 
the moments in Eqs. (| 3 .4411 3 .451) to those in Eqs. (|3.35H3.36l) one finds the equations 

t» = ^—\ Tv +J - K = ^\^\ (3 - 47) 

whose solution yields the value of a = A + iA. It is shown in Ref. [HH that the substitution of Fq(z) given in Eq. 
p.43|) with the value of a satisfying Eqs. (|3.47[) shows that it is indeed the correct solution of the differential equation 
p.39p . Note that, in the context of the BCS theory, the parameters A and A are known as the BCS gap and the 
chemical potential, respectively. Accordingly, Eqs. (|3.47|) are called the BCS gap and chemical potential equations. 

It was already mentioned that the second order moment of the total electrostatic field is related to the energy of 
the corresponding eigenstate. Using Eqs. (|3.18[) . (I3.37[) . (|3.46l) and (I3.47|) . the energy of the eigenstate corresponding 
to the solution given in Eq. (I3.43[) can be found as 

E = - w(2j - k + 1) - Y, 3*o> (l - j^i) + ^ (3-48) 

to the first order in j. Here Er_j\ is given by Eq. p. 111) . 

Note that the neutrino Hamiltonian given in Eq. (|2.23a[) has an overall minus sign with respect to the BCS pairing 
Hamiltonian given in Eq. (12 .26[) . Therefore, although the field Fq(z) given in Eq. (|3.43p yields the ground state of the 
(reduced) BCS pairing Hamiltonian, for the neutrino case, it yields the eigenstate with the highest energy. A study 
of the general properties of Bethe ansatz states in the continuum limit from the perspective of collective neutrino 
oscillations is an open problem and will be considered elsewhere. 



4. THE RANDOM PHASE APPROXIMATION (RPA) 



A large quantum system such as a self interacting neutrino gas can be conveniently studied within an effective one- 
particle approximation whereby it is described in terms of single particles interacting with an average potential created 
by all other particles in the medium. The requirement that the average potential should evolve in line with the time 
evolution of individual particles gives rise to a set of consistency equations. In this section, we formulate the effective 
one-particle description in terms of the RPA method. Our aim is to reach a quantum mechanical effective one-body 
Hamiltonian which describes the interaction of a single neutrino with its background represented by the classical 
polarization vectors. In this picture, the well known time evolution equations of the polarization vectors emerge 
from the RPA consistency condition. In Section 14 A[ we show that the time evolution equations of the polarization 
vectors possess several constants of motion which are basically the expectation values of the many-body constants of 
motion mentioned in Section 12 Bl One of the main results of Section @] involves the diagonalization of the linearized 
RPA Hamiltonian which yields a noninteracting (quasi-particle) basis for neutrinos. We show in Section [4 B I that the 
well known phenomenon of spectral splits can be viewed as a result of the adiabatic evolution of these quasi-particle 
degrees of freedom from flavor to mass eigenstates under appropriate conditions. 



A. Equations of Motion 



In the method of RPA, a quadratic operator 0\Oi is approximated as 

O x 2 ~ O x {0 2 ) + (0!)0 2 - (Oi)<O a ) , 



(4.1) 
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where the expectation values should be calculated with respect to a state \9) which satisfies the condition 

(OiO a > = (Oi)(0 2 ) . (4.2) 

The state \^S) is usually found by solving the resulting RPA consistency equations as illustrated below. 
Application of the RPA method to the neutrino Hamiltonian given in Eq. (|2.23p yields 

H ~ H RPA = ^ujB ■ X + nP ■ J ■ (4.3) 

In Eq. (|4.3p . we introduced the polarization vector which is defined as 

P P = 2<J P > • (4.4) 



The total polarization vectors are defined as in Eq. (|2.4I) . In Eq. (|4.4[) . the expectation values are calculated with 
respect to a state which is assumed to be a tensor product of one-particle states, i.e., 

I*) = l^(Pi)) ® |^(pa)> ® ■ • • ® I^(Pjv)) , (4.5) 

where 

|V(p)>=^e(p)ke>+V'x(p)|^) ■ (4-6) 

The state \ty) defined in Eq. (|4.5|l is an SU(2) coherent state and thus automatically obeys the condition given in 
Eq. (|4.2|) (see Ref. [2l|). Substitution of the state |$) in Eq. (|4.4|) yields the explicit form of the polarization vector 
as follows 4 : 

2Re[rjp)Mp)} 

P p = \ 2Im[r e (p)Mp)} | • (4-7) 

\Mp)\ 2 - \Mp)\ 2 

The linearized Hamiltonian given in Eq. (I4.3[) yields the following Heisenberg equations of motion 

^-X = -i[J u ,H RPA ] = {uB + nP) x t . (4.8) 
ere 

Due to the self consistency requirement of the RPA method, the polarization vectors should then obey 

—P u = (ujB + fiP) x P w (4.9) 

for every u>. Eq. (14.91) can be obtained by taking the expectation values of the isospin operators on both sides of Eq. 
(|4.8|) . Solution of the consistency equations yields the state |\&) postulated at Eq. (|4.5j) . 

Note that in approximating the exact many-body Hamiltonian with the RPA method, one sacrifices the conservation 
of the many-body invariants h u given in Eq. (|2.30l) because the RPA Hamiltonian no longer commutes with these 
operators. However, the expectation values of the operators h u with respect to the state jvP) continue to be invariant. 
The expectation values (h u ) can be easily calculated using Eqs. (|4.2j) and (|4.4[) . The result is 

i u = 2(h u ,) = B-P u + v V IkJjaL. (4.io) 

Z — ' U! — Ul' 



4 An alternative way to define the neutrino polarization vector is to use the density matrix method. For a single neutrino with momentum 
p, the elements of the flavor density matrix are defined as p a fj(p) = (o^(p)oa(p)). For the state given in Eq. 1 14. 51 k this definition 
amounts to 

nM _ ( i<Mp)I 2 r x (p)Mp) \ 
p(p) -\re(p)^(p) iv>*(p)i 2 ) ■ 

Then, one can define the polarization vector though the expansion p(p) = ^ ^1 + 3 ■ P p ^ in terms of the Pauli sigma matrices 3. This 
definition leads to the same result as in Eq. 114,71 1. 
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Here a factor of 2 is introduced in the definition of 1^ for convenience. Using the RPA equations of motion (|4.9[) . it 
is a straightforward calculation to show that 

±J U = 0. (4.11) 

The particular combination 

K a = ^I u = B-P (4.12) 

is already pointed out to be a constant of motion by many authors. Investigation of the N-mode coherence in collective 
neutrino oscillations leads to similar invariants obtained as the scalar product of the polarization vector and its Hilbert 
transform [47]. The invariants we find in the RPA formalism are closely related to the invariants of such N-mode 
coherence. It should be mentioned out that such invariants were first discussed in the context of BCS dynamics [63$ . 

In the mass basis, we denote the polarization vectors by calligraphic letters as in the case of isospin operators, i.e., 
V p = 2(Jp). From Eq. (|2.7|l we see that the conserved quantity given in Eq. (|4.12p is equal to 

K = -V* = . (4 . 13) 

Together with the conservation of (Ni) + (-/V2), this implies that both (N\) and (N2) are individually conserved under 
the RPA evolution. Therefore, although the exact many-body Hamiltonian given in Eq. (|2.23l) strictly conserves the 
number of particles in each mass eigenstate (see Eq. (j2.36[) ). the RPA evolution conserves only their average values. 



B. Stationary Solutions and the Spectral Splits 

It has been pointed out in Section 12 Bl that the BCS pairing Hamiltonian and the collective neutrino oscillation 
Hamiltonian have the same mathematical structure. This equivalence has been exploited in Scction[3]in order to write 
down the exact eigenstates and eigenvalues of the many-body neutrino Hamiltonian. In this section, we will consider 
the eigenstates and eigenvalues of the RPA neutrino Hamiltonian given in Eq. (|4.3I) . These eigenstates correspond to 
the steady state solutions of Eqs. (14.91) . To this end, we will borrow the method of Bogoliubov transformations from 
the RPA-BCS model by which the well known BCS ground state and its excitations are obtained. 

Application of the method of RPA to the BCS pairing Hamiltonian given in Eq. (|2.26p yields 

iJ B R c p * = 5>^-G((T+}T-+T+(T->) . (4.14) 

fc 

In order to find the ground state of the RPA-BCS Hamiltonian for a given number of pairs, one introduces a chemical 
potential through the method of Lagrange variables. In other words, one tries to minimize the Hamiltonian 

- 2 M T* = 2(e fc - UK - G «T+>T" + T+(T~)) . (4.15) 
fc 

The excited states are then obtained from the ground state by applying the excitation operators. Although the 
many-body BCS Hamiltonian _ff B cs conserves the total pair number 2T Z (see Eq. (|2.20|) '). the RPA-BCS Hamiltonian 
given above conserves only its mean value 2(T Z ). Thus, the BCS ground state and the excited states obtained this 
way have a constant mean pair number, but they are not eigenstates of the pair number operator, i.e., they have 
indefinite number of pairs. 

As pointed out at the end of Section [4 Al one has a similar situation for neutrino Hamiltonian and we introduce an 
analogous Lagrange variable which we denote by ui c : 

H RPA + lu c J z = J2("c ~ u)J z + ilP ■ 3 ■ (4.16) 

This way one can find the state which minimizes the RPA neutrino Hamiltonian (|4.3[) for a given value of the constant 
in Eq. (|4. 13|) . Note that, since the Hamiltonian in Eq. (|4.16[) commutes with the total particle number operator N, 
the total neutrino number N is always well defined. Consequently, the Lagrange multiplier in Eq. (|4.16p fixes the 
values of both (Ni) and (Ni). In what follows, we will use a method borrowed from Bogoliubov's solution of the 
RPA-BCS model which yields both the ground and the excited states with these properties [6~j| . 
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Let us begin by pointing out that the Hamiltonian in Eq. (|4. 1 6[) can be diagonalized by rotating the flavor space 
of each energy mode in a different way as described in Section 12 Al Such a rotation can be carried out by the unitary 
operator given in Eq. (|2.8|) . This time we denote the operator by U u and its parameter by z u to emphasize that all 
neutrinos with a given energy are subject to the same transformation while different energy modes undergo different 
transformations. Our aim is to write the Hamiltonian given in Eq. (|4. 16[) in the form 

H RPA + Uc J' = Yl ^UlJ'U K . (4.17) 

UJ 

The first of Eqs. (|2.12j) tells us that for the equality in Eq. (|4.17[) to hold, the rotation parameter z w = e lS " tan# w 
should satisfy 

2A W cos 26> w = (uj c ~uj)+ nV z , 
2\ LO e ±lS " sin2<L = ■ (4-18) 

Solving Eqs. (|4.18[) we find that the phase S u is the same for all energy modes and is given by 

e" = (4.19) 
whereas the rotation angle 9^ depends on the energy and is given by 



'1+ o~\ h sin^ = W- 1 . (4.20) 



V 2 V 2A W J ' V 2 V 2A 

The parameter A w can also be found from Eqs. (|4.18[) as 

A w - \>J(uc - u + vV z ) 2 + v 2 V+V- . (4.21) 

Eq. (|4.17l) motivates the definition of the Bogoliubov transformed quasi-particle operators 5 

ai(p) = C/ ( J ) ai(p)[/„ = cos^ai(p) + e w sin^a 2 (p) (4.22) 
«2(p) = Ula 2 (p)U w = -e~* s sm9 u a%(p) +cosO u a 2 (p) 

together with the corresponding isospin 

J P = UlfpU* , (4.23) 
which we call the a-isospin operators. They are equal to 

J+ = o4(p)a 2 (p) , J p = <4(p)ai(p) , J p = - (a t 1 (p)ai(p) - a 2 (p)a 2 (p)) . (4.24) 

The Hamiltonian given in Eq. (|4.17[) can now be written as the Hamiltonian of a noninteracting system in terms of 
the quasi-particles 

cJ' = 2A - J i = E A - (« f i(p)«i(p) - 4(P)«2(P)) • (4.25) 



In its broadest sense, a Bogoliubov transformation is any mixture of particle operators which preserves the (anti) commutation relations. 
However, the term usually suggests the mixing of creation and annihilation operators because in Bogoliubov's approach to BCS theory 
the electron and hole states are mixed. Although Eq. 114. 2211 is only a rotation in flavor space (apart from the phases) it can also be 
viewed as a special Bogoliubov transformation in the broad sense of the term because the rotated particle operators satisfy the same 
anti-commutation relations as the original ones. In that sense, the distinction between the real particles and the guasi-particles seems 
somewhat arbitrary because one can see the flavor basis as quasi-particles of the mass basis. However the term quasi-particle is usually 
reserved for the noninteracting basis in the presence of an interaction between the particles and the background. Here, our usage of the 
terms Bogoliubov transformation and quasi-particle is due to the established terminology of the BCS theory. 
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The eigenstates of the Hamiltonian given in Eq. (|4.25l) are the ones in which all neutrinos occupy a\ and «2 states. 
In particular, the ground state is an "a 2 condensate." In the language of a-isospins introduced in Eq. (|4.24[) . all 
a-isospins point down in the ground state so that it is equal to the lowest weight state 

\j n ^,-3„^) a = l[4(p) |0) , (4.26) 
p 

where j max = N/2. The states with higher energies can be obtained by exciting the particles from «2 to a\ states 
which amounts to applying the operator J+ defined in Eq. (|4.24l) . In general, the eigenstates and eigenvalues of the 
Hamiltonian in Eq. f|4. 25[) can be written as 

W =\[\iu:,m u ) a and £ = ^2m w A w . (4.27) 

UJ uj 

In Eqs. (|4.26l) and (I4.27p . the subscript a points to the fact that these are the eigenstates of the a-isospin 6 . The 
ground state given in Eq. (|4.26|) is the special case of Eq. (|4.27j) in which — ~j u for all energy modes. 
For consistency, the states in Eq. (|4.27p should satisfy 

2{^\J+\i>) =V+ and 2{tl}\J'\il>) =V Z . (4.28) 

It is an easy exercise in algebra to show that 

2^\J + |V.) = 4 J2 ™u> — ^ = e a n\V + \ £ ^ (4.29) 
« 1 + |^| 2 V A - 

and 

In this case, the consistency equations (|4.28|l become 

2pL ~ ^ 2A W 2 M ~ 2s 2A U ' 1 ' ' 

UJ UJ 

The consistency equations given in Eqs. (|4.9[) differ from those given above in that the former ones are valid for a 
generic time dependent state whereas the later ones are valid for the eigenstates (i.e., stationary states) of the RPA 
Hamiltonian. At this point, we would like to note the similarity between Eqs. (|3.47[) and Eqs. (14.31[) . From Eq. (|4.27[) 
we see that the eigenstate with the highest energy has m w = j u . For this particular state, the consistency equations 
(|4.31[) are the same as those given in Eqs. p. 471) with the substitution A = lo c — f/P z and A = ^|"P + | together with 
K = j + V z /2 (see Eq. (JS^gj). 

A special but very interesting case of Eqs. (|4.31[) is obtained in Ref. [2(| in connection with the spectral splits of 
neutrinos which occur if the parameter \i defined in Eq. (|2.22j) changes in an adiabatic way from oo to (also see 



It is evident from Eq. 1 14. 23H that these states are related to the corresponding eigenstates \ju ,m UJ ) m of the mass isospin by 

|jw)f^o))a — U^jwimuj)™ t 
which in turn are related to the eigenstates |ju,m w )j> of the flavor isospin by 

\ju,mu)m = U<\j u ,mu}f ■ 
Therefore, we can write the eigenstates in Eq. 114. 2711 as 

W = EI b"> m ")i = U lYl \0u,,m u )m = UlU ] Y[ \0u»m u )f ■ 

UJ UJ UJ 

Note that the transformation described by the operator U w is not a global rotation. From Eq. 114. 201 1, it is clear that each energy mode 
ui undergoes a different rotation. As a result, this transformation does not preserve all scalar products. In particular, the total o-isospin 
quantum number of a state is not always the same as its total mass or flavor isospin quantum number. The state in Eq. 14.2611 . for 
example, does not live in the j max representation of mass or flavor isospin although it lives in the jr max representation of the a-isospin. 
However, j u is a well defined quantum number in all three bases because, within each energy mode ui, the transformation induced by 
U^> is a "global" one. 
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Ref. [23[ for an in-depth review). In this case, the Lagrange multiplier lj c introduced above plays the role of the split 
frequency. Here we will show that this phenomenon can also be understood in terms of the adiabatic evolution of 
quasi-particle states. Namely, we will show that when the neutrino density is high, the quasiparticle states coincide 
with flavor eigenstates whereas in vacuum they coincide with mass eigenstates. As a result, neutrinos emerging in 
flavor eigenstates in the high density region gradually transform into mass eigenstates if the neutrino density decreases 
in an adiabatic way. The Lagrange multiplier uj c determines the final energy distribution of mass eigenstates. 

In order to see how this phenomenon takes place, let us assume that initially the neutrino density is so high that 
the — > oo limit is practically attained and that all neutrinos initially occupy flavor eigenstates v e and v x . Note that 
here it is not assumed that all neutrinos occupy the same flavor state. Neutrinos emanating from the surface of a 
proto-neutron star after a supernova explosion are believed to meet these initial state assumptions. In /i — > oo limit, 
Eq. (I4.21[) gives 2A W = fj\V\ — fi\P\ for all neutrinos. Since we assumed that all neutrinos are in flavor eigenstates, 
we have P+ = P~ = in which case Eqs. ([2~T4]) lead to V z = cos26P z . Substituting these in Eq. (l4T2"0")) we find 



cos6L = 



\ 



1 pz 

- 1 + — cos 261 

2 \ \P\ 



\ 



1 / P z 

- 1 - — cos 29 

2 I IPI 



(4.32) 



Now, let us assume for a moment that we have more u e than v x so that P z > 0. In this case, Eq. (|4.32p yields 
cos 6^ = cos 9, i.e., the angle 9 U is equal to the vacuum mixing angle for all neutrinos. Also note that Eqs. f!2 . 14[) lead 
to V + = sm29P z which tells us that the phase S defined in Eq. (I4.19|) is equal to zero for this particular case. These 
observations lead to the result that z u — z = tan0 for all neutrinos. In this case, the operator simply becomes 
equal to the operator U defined in Eq. (|2.11l) . As a result, one arrives at 



ai(p) = U^ai_(v)U = a e (p) and a 2 (p) = L/ t a 2 (p)C/ = a x (p) , 



(4.33) 



where we used Eq. (|2.10l) . The result stated in Eq. (|4.33p is important because it tells us that in the limit where 
fi — > oo, the quasi-particle operators become equal to the flavor particle operators and thus the state in which all 
neutrinos occupy flavor states is an eigenstate of the RPA Hamiltonian (I4.25[) in this limit. Note that, although we 
assumed that we have initially more v e than v x , the opposite scenario can be worked out in a similar way and shown 
to lead to the same conclusion with ai(p) = — a x (p) and a 2 (p) = a e (p). 

Since one of the states in Eq. (|4.27l) represents the initial state of the system, then we have a steady-state case 
with no time evolution as long as fj, stays constant. It is already known that in the limit where /x — > oo the state in 
which all neutrinos occupy flavor eigenstates is a steady-state solution which is consistent with the above discussion. 
However, if /i decreases (which happens when the neutrinos emanate from a source), then the quasi-particle operators 
ai(p) and Q!2(p) change accordingly and the initial state cannot remain as an eigenstate. However, if the change in 
fi occurs in an adiabatic way, then the state evolves steadily in such a way that it continues to keep its original form 
in terms of the quasi-particles while the quasi-particle operators «i(p) and a2(p) themselves slowly change. In the 
limit where \x — > 0, Eq. (|4.21[) yields 



2\ui = 



and Eqs. P~2"0"|) lead to 




sin 9 U 



'1 



f 



UJ 



(4.34) 



(4.35) 



\uj c - u\j " y 2 

As a result, in \x — > limit one obtains cos 9^ — 1 and sin 9^—0 for those neutrinos with lj < uj c leading to 

ai(p) = ai(p) , «2(p) = a 2 (p) ■ (4.36) 

On the other hand for those neutrinos with lj > lo c1 the /x — > limit yields cos6* w ~ and sinf^ = 1 leading to 

«i(p) = «2(p) , a 2 (p) = ai(p) , (4.37) 

up to some inconsequential phase factors. We see that, a state which starts off with all neutrinos in flavor eigenstates 
in the [i — > oo limit can adiabatically evolve towards the /i — > limit in such a way that those neutrinos with uj < uj c 
and those with uj > uj c end up in the opposite mass eigenstates. For example, if we initially have more v e than u x , 
then those electron neutrinos with lo < u c evolve into the first mass eigenstate while those with uj > oj c evolve into 
the second mass eigenstate. In this case, the opposite is true for v x . Such an evolution leads to a spectral split as 
discussed in Ref. 1201 - 
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In order to find the split frequency u> c , one should solve the consistency equations given in Eqs. (I4.3ip . Since we 
assumed that initially all neutrinos are in flavor states, one has Im^ = P* = e w 17-^1 where e u = ±1 depending on 
whether P* initially points up (more v e ) or down (more v x ). As a result, the consistency equations given in Eq. (|4.31[) 
reduce to 

1 = T , ^ and u, c = Y ; (4.38) 

These equations were also obtained in Ref. [20| where a method based on the time evolution of polarization vectors 
is used. The derivation presented here in terms of the evolution of quasi-particles is physically equivalent to the one 
given in Ref. [H. 



5. ANTINEUTRINOS 



In this section, we include antincutrinos in the formalism introduced in previous sections. Let us denote the fermion 
operators for a neutrino with momentum p in the flavor state v a by a a (p) and in the mass state Vi by a, (p) where 
a = e, x and i = 1, 2. It has been pointed out by many authors that describing antineutrinos with the rotated spinor 
(— instead of the regular spinor (p\,i?2) is more advantageous (see, for example, Refs. 0, Hl|) because this 
leads to a seamless integration of antineutrinos into the formalism making the physics more transparent (see Appendix 
iBj) . In order to implement this, one can first define the antineutrino isospin operators as in Eqs. (|2.2p and (|2.7p and 
then perform a transformation as in Eqs. (|2.9I) with 9 p = ir/2 and <5 P = for all momentum modes 7 . This leads to 

ai(p) = -a 2 (p) , 5 2 (p) = oi(p) (5.1) 

and similarly for the particle operators of flavor states. We denote the corresponding flavor and mass antineutrino 

isospin operators by J p and J v which are defined as in Eqs. (I2.2[) and (I2.7[) . respectively. 

With the introduction of antineutrinos, we modify the summation convention introduced in Eq. (12.41) as 



A,= ^A p , A_ H = and A = J2 A " ( 5 - 2 ) 



|p|=p IpI=p 

and also introduce 



5m 



(5.3) 



for antineutrinos. Note that the last sum in Eq. (15.21) over uj now includes the negative values for antineutrinos as 
well as the positive values for neutrinos. 

The Hamiltonian describing the vacuum oscillations and self interactions of neutrinos and antineutrinos in the single 
angle approximation is given by Eq. (|2.23l) , except that now the new convention introduced in Eq. (|5.2p applies [56| . 
As a result, all of the main results of this paper, including the invariants given in Eqs. (|2.30p and (|4.10p . can be 
generalized to include antineutrinos by extending the sums over uj to include negative values. We only would like to 
make a few comments in what follows. 

With the inclusion of antineutrinos, the number of invariants is now doubled because one has the invariants h u for 
neutrinos and the invariants h— u for antineutrinos. Considering the diagonalization of the Hamiltonian described in 
Section (13 Ap . j now represents the total isospin quantum number of both neutrinos and antineutrinos in accordance 
with Eq. (|5.2p . In the case of the electrostatic analogy, one has additional negative charges — j_ u located on the 
negative x-axis at positions — uj. In this case, one can still apply the formalism described in Section |3 Bl by shifting 



7 Note that for 8 P = it/2 the parameter z p is singular in Eq. 112.81 1. However, the transformation is still well defined because the operator 
Up can be written in several alternative forms. One of them, which is not singular for 8 P = 7r/2, is given by 

JT - „E„ ^(J* sin S p +jy cosip) 
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the origin of the coordinate system in Fig. [TJ Finally we note that when antineutrinos are treated as described above, 
their polarization vectors defined as in Eq. (|4.4I) are given by 



This definition involves an overall minus sign in comparison to the one which is frequently used in the literature. 



In this paper, we studied the symmetries and the associated constants of motion of the collective neutrino oscillation 
Hamiltonian by taking into account both the vacuum oscillations and the self interactions of neutrinos. We examined 
the system both from the exact many-body perspective and from the point of view of an effective one-body description 
formulated with the application of the RPA method. We showed that, under the single angle approximation, both the 
many-body and the RPA pictures possess many constants of motion manifesting the existence of associated dynamical 
symmetries in the system. 

The existence of these constants of motion make the system completely integrable in the sense that the exact 
eigenstates and eigenvalues can be found in an analytical way. We wrote down these eigenstates and eigenvalues for 
both the exact many-body and the RPA Hamiltonians. In the case of the many-body Hamiltonian, the eigenvalues 
and eigenstates were found with the application of the Bethe ansatz method and they depend on the solutions of 
Bethe ansatz equations which are analogous to the equilibrium conditions of an electrostatic system in two dimensions. 
In the case of the RPA Hamiltonian, the eigenstates and eigenvalues were found by applying a suitable Bogoliubov 
transformation which brings the Hamiltonian into the form of a noninteracting system in terms of quasi-particles. We 
wrote down the consistency equations for all RPA eigenstates and showed that these equations reduce to a particular 
but physically very important case which was studied earlier in Ref. [201 ] in connection with the spectral splits of 
neutrinos. We showed that the spectral splits can be understood as the adiabatic evolution of quasi-particle states 
from a high density region where they coincide with flavor eigenstates into the vacuum where they coincide with mass 
eigenstates. 

In general, the existence of constants of motion offers practical ways of extracting information even from exceedingly 
complex systems. Even when the symmetries which guarantee their existence is broken, they usually provide a 
convenient set of variables which behave in a relatively simple manner depending on how drastic the symmetry 
breaking factor is. In this paper, we omitted an ordinary matter background and concentrated on a neutrino gas 
which undergoes vacuum oscillations as well as self interactions. Although an ordinary background of p, n, e~ and 
e + coexists with neutrinos in most astrophysical sites, it is only the net electron fraction which plays a part in flavor 
evolution of neutrinos through a diagonal potential in flavor basis which is proportional to y/2G F\n e - —n e + \ . Neglecting 
this potential can only be justified when the asymmetry of the electron background is small in comparison to the 
neutrino background. On the other hand, even under the circumstances where the electron background asymmetry 
cannot be ignored, the constants of motion presented in this paper can still be useful as a set of convenient variables 
because, since they commute with the neutrino gas part of the Hamiltonian, their time variation only come from their 
commutator with the electron background term. 

As was recently illustrated in Ref. [13] , existence of such invariants naturally leads to associated iV-mode collective 
neutrino oscillations. These collective oscillations are closely related to the m-spin solutions presented in Ref. [§3[ in 
the context of the BCS model. However, symmetries alone do not guarantee the stability of such collective behavior. 
Some issues related to the stability of the m-spin solutions of the BCS model have been addressed in Refs. [63l. l65l|. 
In the case of the collective neutrino phenomena associated with our invariants the question of stability is an open 
problem and could be illuminated by numerical studies. 
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Appendix A: Alternative Ways of Expressing Quantum Invariants 

In this appendix, we will briefly mention two alternative ways to rewrite the constants of motion mentioned in the 
text. The first method is based on the following sum of the many-body invariants h u and L u : 

C n + ^iw"- 1 ^ . (A.l) 

UJ UJ 

Clearly C„ is a constant of motion for every value of n but here we specifically assume that n is a positive integer. If 
we define the quantities 

Q n = 5>"jL, (A.2) 

UJ 

then C n defined in Eq. (|A.1[) can be written as 

C n = B ■ Q n + yt (Qn-l ■ <2o + Qn-2 ■ Ql + ■ ■ ■ + Ql ■ Qn-2 + Qo • Qn-l) ■ (A.3) 



The invariant Cq was already mentioned in Eq. (|2.32[) . The Hamiltonian itself is equal to C\ as can be seen from Eq. 

A similar method can be used to rewrite the constants of motion of the RPA formalism. One starts from the sum 
of the RPA invariants I u and |P W | 2 given by 

which is clearly an invariant for every n. Defining the quantities 

P n = £>"P w , (A.5) 

UJ 

it can be shown that K n is equal to 

K n = B ■ P n + | (P„_! • P + Pn-2 ■ P\ + • " • + Pi ■ Pn-2 + P • P„_l) . (A.6) 

Another alternative way to express the many-body quantum invariants, which is frequently encountered in the 
literature on Gaudin formalism, is based on the following sum: 

*w~EsVE&V (A - 7) 

UJ UJ v ' 

Here A is a complex valued parameter and H (A) is clearly a constant of motion for every value of A. Defining 

2(A)- E^T (A-8) 
^— ' ui — A 

UJ 

and using Eq. (|2.30l) . one can show that the operator H(X) defined in Eq. (|A.7[) can be written as 

ff(A) = Q z {\) + mQ(A) • Q(A) . (A.9) 
H(X) is said to form a one parameter family of conserved quantities. It satisfies 

[H(X),H(n)]=0 and [H(X),H]=0 (A.10) 

for every complex value of A and \x. 
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In the case of the RPA formalism, one similarly defines the sum of the RPA invariants I u and |fL| 2 as follows: 
-D(A) is a constant of motion for every complex value of A. Defining the vectors 



cj — A 

one can show that it can be written as 



£>(A)=P*(A) + |p(A)-P(A) • (A.13) 
D(A) also forms a one parameter family of conserved quantities. 

Appendix B: Antineutrino Hamiltonian 

The effective many-body Hamiltonian describing the collective flavor evolution of neutrinos in a dense medium in 
the presence of neutrino self interactions was obtained in Ref. [56] by expanding the full neutral current interaction 
and keeping only the forward scattering terms. The underlying group is SU(3) x SU(3) representing all three neutrino 
and antineutrino flavors. Reducing to the SU(2) x SU(2) subgroup for the two neutrino case yields the following 
effective many-body Hamiltonian: 

H = u I - cos2g J° + sin2fl J " + J " _ C os2fl J\ + sin 20 + ^— ) (B.l) 



ui>0 



^ ) (1 cos ^pq) ' Jq + Jp " Jq 2J p J q Jp Jq Jp Jq 



Here J are the flavor isospin operators for antineutrinos defined as in Eq. (|2.2|) . They are subject to the same 
summation convention as A given in Eq. (|5.2p . The terms involving the neutrinos and antineutrinos in the Hamiltonian 
(IB.ip are dissimilar. However, if one performs the transformation described at the beginning of Section [5] and defines 
the antiparticle operators a as in Eq. (I5.1[) . then the corresponding isospin operators are related by 

(Jp, Jp", Jp) = (~Jp, ~Jp, -Jp") ■ (B-2) 
In terms of the transformed antineutrino isospin operators, the Hamiltonian in Eq. (jB.ll) becomes 

h = cos26> J ° + sin2d J " + J " + cos26 ^ ~ sin29 J ~" + I ^ B - 3 ) 

+ — cos ^pq) (jp ■ jq + Jp ■ Jq + 2jp ■ j 

If we define the oscillation frequencies of antineutrinos with a minus sign as given in Eq. (|5.3[) . then (in the single 
angle approximation) the total Hamiltonian (IB.3|) takes the form of Eq. (|2 . 23[) with the sum over uj extended to 
include the negative values. 
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